Dynamical Cluster Quantum Monte Carlo Study of the Single Particle Spectra of 

Strongly Interacting Fermion Gases 
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We study the single-particle spectral function of resonantly-interacting fermions in the unitary 
regime, as described by the three-dimensional attractive Hubbard model in the dilute limit. Our ap- 
proach, based on the Dynamical Cluster Approximation and the Maximum Entropy Method, shows 
the emergence of a gap with decreasing temperature, as reported in recent cold-atom photoemission 
experiments, for coupling values that span the BEC-BCS crossover. By comparing the behavior of 
the spectral function to that of the imaginary time dynamical pairing susceptibility, we attribute 
the development of the gap to the formation of local bound atom pairs. 



One of the most exciting recent developments in cor- 
related systems has been the observation of pairing and 
superfluidity of ultracold atomic fermions interacting via 
an s-wave Feshbach resonance [IHH . The precise tunabil- 
ity of the Feshbach resonance mediated interaction allows 
the exploration of superfluidity for coupling extending 
from the BCS regime, with loosely-bound Cooper pairs, 
to the BEC regime, with local pairs at strong attraction. 
Thus such experiments provide a highly controllable set- 
ting to study phenomena that appears in a range of other 
systems, including neutron stars quark matter 
and correlated electronic systems such as the high tem- 
perature superconductors [6|. 

While the low-temperature properties of ultracold 
fermion gases are relatively well understood and essen- 
tially described by the BEC-BCS crossover wavefunc- 
tion [i| , the finite temperature behavior still presents sev- 
eral mysteries, particularly in the unitary regime, where 
the s-wave scattering length a s between the two species 
of fermion diverges. At strong couplings preformed pairs 
are predicted to exist above the superfluid transition tem- 
perature y| , possibly yielding a pseudogap in the density 
of states 7|. But, how the pseudogap phenomena is re- 
flected in experimental probes, such as radio frequency 
spectroscopy [8|, |9( , is still an open question. 

To shed light on the spectral properties of resonantly- 
interacting fermion gases, we compute the single particle 
spectral function A(k, uj) as a function of frequency u) and 
wavevector k. The spectral function has recently been 
measured in experiments with trapped 40 K atoms via a 
novel momentum-resolved generalization fiol j of radio fre- 
quency spectroscopy @, Q . Various analytical methods 
have been developed to compute A(k, uj) in the BEC-BCS 
crossover based on the ladder approximation ll|, Il2j . 
self-consistent many-body theories [pj - [l5| and large-A^ 
expansions [l6|. These approaches are hampered by the 
lack of a small parameter in the unitary regime. 

In this letter we combine the Dynamical Cluster 
Approximation 17, [l8| and the Maximum Entropy 
Method [l9| to obtain the spectral function in the ther- 
modynamic limit without uncontrollable approximations 
(see also Ref. [13 )■ We observe a strong suppression of 



low-energy spectral weight in conjunction with the ex- 
pected onset of pairing correlations. Our simulations 
also capture the concomitant formation of a double-peak 
structure in ^4(k, uj) with decreasing temperature in the 
strong-coupling regime. At unitarity, we find a momen- 
tum distribution «/,. consistent with the predictions of 
Tan [211 ] . and an estimate of the universal parameters £ 
and £ characterizing the unitary Fermi gas [Uj]. 

The Dynamical Cluster Approximation (DC A), unlike 
conventional QMC, allows the extraction of the single- 
particle spectral function A(k, uj) in the thermodynamic 
limit. One disadvantage of the DCA in the present set- 
ting is that it applies to lattice models, while cold atoms 
possess the itinerant-fermion dispersion et = k 2 /2m. 
However, we address this by focusing on the dilute limit 
of the three-dimensional cubic Hubbard model 
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where c,- CT annihilates a spin-cr fermion at site i, the angle 
brackets restrict the sum over nearest neighbor sites with 
the matrix element t, and U is the attractive on-site in- 
teraction. Here, the chemical potential fi is tuned to keep 
the filling at the small value (n) = 0.3, where (n) = 1 is 
half filling. The Fermi energy associated with this filling 
is Ep — 3.414i. In this dilute limit, the hopping dis- 
persion £k = 2f [3 — J2i= x y z cos k»] can ^ e approximated 
by its low k expansion ej, ~ tk 2 , thus describing free 
fermions if we take t — 1 /2m. 

By tuning U close to U c — — 7.92i, the interaction 
strength where a two-body bound state forms in vac- 
uum [23| , we can capture the unitary regime in cold Fermi 
gases. We expect that interaction effects are insensitive 
to the underlying lattice when the scattering length 

as - 87Tt(\u\ - \u c \y {Z) 

exceeds the lattice spacing. We find, as discussed below, 
that Eq. (JlJ works well in the BCS and unitary regime, 
but begins to deviate from the free-fermion results in 
the deep-BEC regime of tightly-bound local pairs where 
lattice effects cannot be neglected. 
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FIG. 1: (color online) The spectral function A(k,uu) of the three-dimensional attractive Hubbard model with filling (n) = 0.3, 
at temperature T/Ef = 0.049, and several coupling strengths, (a) On the BCS side of the Feshbach resonance, U = — 6.88i, 
l/(a s fcp) = —0.99. (b) At unitary coupling, U = U c — — 7.92i, l/(a s kp) — > 0. (c) On the BEC side of the Feshbach resonance, 
U = -8.4t, l/(a B fcp) = 0.38. (d) In the BEC region, (/ = -9.2t, l/(a s k F ) = 0.92. In each panel, the white curve is the 
noninteracting dispersion ek — /U- 



The DCA covers the lattice with clusters of iV c = I? c 
sites using periodic boundary conditions. Correlations 
up to length L c are considered exactly while long range 
correlations are treated at the mean-field level [lj], ■ 
We use the Hirsch-Fye quantum Monte Carlo (QMC) 
algorithm [24| as the quantum solver 25 1, and the max- 
imum entropy method [l9| to analytically continue the 
imaginary-time spectra to real frequencies. In the DCA 
the momentum resolution of the self energy £ is deter- 
mined by a set of N c i^-points corresponding to the clus- 
ter reciprocal space 18| . We use a Fourier Transform in- 
terpolation scheme to approximate E on other k-points 
in the Brillouin zone. We use a Betts cluster [26[ with 
N c = 24. Our Hirsch-Fye QMC algorithm is efficient as 
long as \U\ is below the bandwidth 12i, allowing us to 
study the BEC-BCS crossover but not the asymptotic 
BEC limit. 

We have obtained the single-particle spectral func- 
tion A(k, v) = ImG(k, w), with G(k, uj) the retarded 

Green's function, across the BEC-BCS crossover and for 
a wide range of temperatures. Our displayed results are 
along the line connecting the center of the first Brillouin 
zone r = (0,0,0) to the corner K = (7r, tt, 7t), with 
k referring to the wavevector magnitude. The spectral 
function along other directions is qualitatively equivalent. 
We focus on the large scattering length limit with — 1 < 

< 1, where we approximate the Fermi wavelength 

k F a s 

k-p by the itinerant-fermion formula kp = {"S^n) 1 ^ . Our 
simulations are in the normal state, meaning we cannot 
fully describe the broken symmetry state below the su- 
pcrfluid transition at T c . Above T c , however, expect our 
results to be quantitatively accurate. Previous work 23] 
has demonstrated that Eq. ([IJ yields reliable results for 
the transition temperature of itinerant superfluid Fermi 
gases. By taking the extrapolated transition tempera- 
ture T c from Ref. we get T c ~ 0.064^ at unitarity 
coupling and filling (n) = 0.3. 



The most striking feature of our results, shown in 
Fig. [U is the strong suppression of the spectral weight 
near zero frequency with increasing attraction towards 
the BEC limit. Fig. Q] displays results for four cou- 
pling values at the temperature T/Ep = 0.049. In all 
plots, the behavior at high-fc is an artifact of the un- 
derlying lattice that we shall ignore; this is reasonable 
since experiments only probe the occupied part of the 
spectral function. Fig. [TJi, in the weak-coupling BCS 

regime ( = —.99 ) , shows an essentially free-particle 

\kpa s J 
dispersion for k < k-p and u> < 0. The narrowing of 

the linewidth as the Fermi surface is approached from 
k — » fcp is clearly seen, consistent with the system being 
a Fermi liquid in this regime. With increasing coupling a 
two branch structure begins to form at the Fermi surface. 
While this behavior is barely visible at the unitary point, 

= 0, shown in Fig. [T](b), it increases in magnitude 

with increasing interaction strength, as seen in panels 
(c) and (d). In Fig-H^d), the deep BEC limit, the two 
branches of the spectral function are clearly visible, and 
qualitatively consistent with the recent cold-atom pho- 
toemission experiments flpj ] . 

To examine the onset of the pseudogap with decreasing 
temperature, we plot A(k,u>) at fixed k as a function of 
energy in Fig. [2] For each coupling we have chosen the 
value of momentum k for which the distance between the 
two spectral branches is the shortest. We find the pseudo- 
gap, with decreasing T, only for coupling values beyond 
the unitary point. In the BEC regime, as seen in Fig. [5J 
the high temperature spectral functions show a Gaus- 
sian shape at high temperatures that develops a pseu- 
dogap at lower temperatures. The pseudogap separates 
two asymmetric peaks, reminiscent of recent analytic re- 
sults [l2(. We define the pseudogap A as one-half the 
minimum peak-to-peak distance between the two spec- 
tral function peaks; our values for A at the lowest tem- 
perature (characterizing the evolution of the pseudogap 
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FIG. 2: (color online) Plots of A(k, cj) at fixed momenta as a 
function of frequency at successive decreasing temperatures, 

for a) Close to unitarity ( = 0.38 ] and b) in the deep 

V k F a a J 

BEC regime ( -: = 0.92 ] . Since we expect superfluidity 

\k F a s J 

only below T c < 0.064_Ef [23|| . the suppression of spectral 
weight above this temperature represents the pseudogap. 
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FIG. 3: (color online) (a) Chemical potential /x as a function 
of coupling at T/Ef = 0.049. (b) The pseudogap A extracted 
from A{k,ui) (square symbols), (c) Momentum distribution 
n(k) for different values of l/(a s fcp). The dashed lines are 
used to guide the eye, and the unitarity curve has been fit at 
large k to the Tan relation [2lj . 

whether the Tan relation [2l| 

n(k) = C/k 4 , for k -> oo, (5) 



across the BEC-BCS crossover) are plotted in Fig. |3|b). 

Even though our simulation does not explore the su- 
perfluid state, it is able to capture the onset of pairing 
correlations and estimate the universal parameters £ and 
£ characterizing the unitary Fermi gas |22j. For example, 
the internal energy per particle can be expressed as 



E/N = -E F [£ 



c 



Oi 



1 



(k F a s 



r)]- 



(3) 



near unitarity |22l . |27[. We quantify £ from our ex- 
tracted values of [i, that is adjusted to always keep filling 
(n) = 0.3. Fig. [3ja) displays fi as a function of coupling. 
The chemical potential is positive and close to Ep in the 
BCS regime and reduces monotonically with increasing 
attraction towards the BEC limit of local pairing. We 
find fi/Ep = 0.31 at unitarity, yielding our prediction 
for the universal parameter n/Ep = £, = 0.31, which is 
consistent with recent experimental and theoretical re- 
sults, see e.g. Refs. dill]. 

We also calculate the momentum distribution 



n(k) 



djjjrvp{oS)A{k, w), 



(4) 



1 



the Fermi function. The momen- 



with n-p(u) — 

turn distribution has been previously measured in the 
BEC-BCS crossover [2{|. Fig. [3jc) displays our results 
for n(k) showing the expected behavior: a sharp step in 
n(k) in the BCS regime that broadens with increasing 
attraction towards the BEC limit. To gain additional in- 
sight, we fit these curves to cq + c\/(\k/k-p\ + C2) 4 to test 



governing the large k behavior of the momentum occupa- 
tion, holds for our system. Although our system mimics a 
free-fermion dispersion at small momenta, at large k lat- 
tice effects become inescapable, and it is not clear if the 
Tan relation will be valid since it probes the extremely 
short distance behavior of interacting fermions. Despite 
this, we find that, at unitarity, the optimal fitting pa- 
rameters are cq = 0.019, c\ = 0.148 and c-i = —0.057. 
The Co and C2 are reasonably small and the fit, shown as 
a solid curve for the unitary plot in Fig. [3tc), is seen to 
be quite accurate, showing that our system does exhibit 
the universal properties of a unitary Fermi gas. From the 
theoretical relation c\ = 2^/5n, with £ governing the in- 
teraction dependence of the energy per particle near uni- 
tary [2l| (Eq. Q), we extract the prediction £ ~ 1.16. 
In the BCS and BEC regimes, the best-fit values of cq 
and C2 are larger, so that we do not find a good fit to the 
Tan formula. In particular in the deep BEC region, one 
expects significant occupation of high momentum states. 
While this certainly occurs in our system, the finite size 
of our momentum-space Brillouin zone restricts the possi- 
ble occupied momenta, preventing our system from truly 
accessing the deep BEC limit. 

It is natural to assume that the pseudogap observed 
in our single-particle spectra arises from the formation 
of s-wave pairing correlations. To test this assumption, 
we study the local s-wave pair susceptibility X s (i?, t) at 
spatial distance R = i — i' and imaginary time r. Setting 
R = and r = (3/2, we have 

X s (0,/3/2) = Ij;(Al(W/2)A,(t' 1 0)« ! , j ,+li.c.) 1 (6) 
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FIG. 4: (color online) Local on-site pairing susceptibility 
X"(-R = 0, /8/2) (a) as function of temperature T for different 
couplings and (b) as function of coupling U for T/Ef = 0.049. 

where A s (i,r) — Ci^(r)cii(T) is the s-wave pair opera- 
tor and the average occurs over system sites. This quan- 
tity, X s (0> /?/2), probes local pairing correlations with the 
longest distance along the imaginary time direction. 

In Fig. @J we plot x s (®i as function of temperature 
(left panel) and coupling (right panel). On panel (a), we 
can see that x s (0,/3/2) vanishes when T — > at weak 
coupling. This signals that, in the weak coupling region, 
atoms can not form stable local pairs; although, of course, 
Cooper pairing correlations will set in at very much lower 
temperatures, outside the range of our computational re- 
sults. When \U\ is large enough, towards the BEC limit, 
X s (0,/3/2) shows a finite tZ-dependent intercept in the 
zero temperature limit, showing the establishment of lo- 
cal pairs. Fig. H£b) shows the coupling dependence 
of x s (0, /3/2) at low temperatures. While this quantity 
is small at weak coupling in the BCS regime, near uni- 
tarity it undergoes a rapid increase as local pairing is 
established. 

To summarize, we have demonstrated that the Dy- 
namical Cluster Quantum Monte Carlo, along with 
the Maximum Entropy Method, may be used to make 
quantitatively-reliable predictions for the single-particle 
spectral function of unitary Fermi gases. Our spectral 
function results capture features seen in recent 40 K pho- 
toemission experiments [lpj . We find the onset of a 
pairing pseudogap with decreasing temperature (above 
the expected superfluid transition temperature) for cou- 
pling values beyond the unitary point in the BEC regime. 
Since our method treats the short range correlations of 
the system exactly, it holds particular promise for ob- 
taining quantitatively accurate predictions for radio fre- 
quency spectra and photoemission data in the strongly- 
interacting unitary regime of superfluid Fermi gases. Fu- 
ture work will extend this method to incorporate the ef- 
fect of a background trapping potential, a necessary step 



to quantitatively understand experimental results. 
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